Mon. Not. R. Astron. Soc. 000. [THTT1 T2013) Printed 18 April 2013 (MN Wl&i style file v2.2) 



Constraints on Warm Dark Matter models from 
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ABSTRACT 

Structures in Warm Dark Matter (WDM) models are exponentially suppressed below 
a certain scale, characterized by the dark matter particle mass, m x . Since structures 
form hierarchically, the presence of collapsed objects at high-rcdshifts can set strong 
lower limits on m x . We place robust constraints on m x using recent results from 
the Swift database of high-redshift gamma-ray bursts (GRBs). We parameterize the 
rcdshift evolution of the ratio between the cosmic GRB rate and star formation rate 
(SFR) as oc (l+z) Q , thereby allowing astrophysical uncertainties to partially mimic the 
cosmological suppression of structures in WDM models. Using a maximum likelihood 
estimator on two different z > 4 GRB subsamplcs (including two bursts at z > 8), we 
constrain m x > 1.6-1.8 keV at 95% CL, when marginalized over a flat prior in a. We 
further estimate that 5 years of a SVOM-like mission would tighten these constraints 
to m x > 2.3 keV. Our results show that GRBs are a powerful probe of high-redshift 
structures, providing robust and competitive constraints on m x . 
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1 INTRODUCTION 

The current concordance cosmology, in which structure for- 
mation proceeds in a hierarchal manner driven by pres- 
sureless cold dark matter (CDM), has been remarkably 
successful in explaining the observed properties of large- 
scale structu res in the Universe (e.g.. iTegmark et al.ll200a : 
iBensqrfeOld) and the cosm ic microwave background (CMB) 
(e.g.. Komatsu et al.ll201~l] ). Such observables probe scales 
in the range ~ 1 Gpc down to ~ 10 Mpc. On smaller 
scales, < 1 Mpc, there are still some discr epancies be- 
tween standard ACDM and observations (e.g., iMenci et all 
l2012t ). For instance, N-body simulations predict more satel- 
lite galaxies than are observed both ar ound our g alaxy 
(the s o-called "miss ing s atellit e problem" : iMoore et all e.g.. 
1 19991 : iKlvpin et all e.g., Il999h . and in the field as recently 
noted by the ALFA LFA survey (e.g.. |Papastergis et al.ll201ll : 
iFerrero et aDl2012h . Furthermore, simulations of the most 
massive Galactic CDM subhaloes are too centrally con- 
densed to be consistent with t he kinematic data of th e brigh t 
Milky Way satellites (e.g., iBovlan-Kolchin et all l201ll ). 
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Moreover, observations of small galaxies show that their 
central density profile is s hallower than predicted by CDM 
N-body simulations (e.g. 



Moore 19941: [de Blok et al 



2001 



Donato et ail 120091 : [m accio et all 120121 : iGovernato et all 
2012). 



Baryonic feedback is a popular prescription for re- 
solving such discrepancies. Feedback caused by supernovae 
(SNe) explosions and heating due to the UV background 
may suppress the baryonic content of low-mass haloes (e.g., 
Governato et al . 2007; Mashc henko et al . 2008; Busha et all 
20101 ; ISobacchi fc Mesingerl l2013bj). and make their in- 
ner density profile sh allower (e.g.. Ide Souza fc IshidalfeoiOl ; 
Ide Souza et allfeoill ). However, accurately matching obser- 
vati ons is still difficult even wh en tuning feedback recipes 
(e.g. IBovlan-Kolchin et allfeoi^ . 

An alternative explanation might be found if dark 
matter (DM) consisted of lower mass (~ keV) parti- 



cles and thus was "warm" ( WDM; e.g., Bode et all 
20011; iKhlopoy fc Kouvarisll2008l; ^ 



de Vega et al 



Kamada et all 



2012 



Kang et all 



de Vega fc Sanchezl 
20131 ; bestri et all 



201S 



2015 



2013). The resulting effective pressure and 



free-streaming would decrease structure on small-scales, 
though again fine-tuning might be required to fully 
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match all the obse r vations (e.g. Bovlan-Kolchin et al.|[201ll ; 
iMaccio et al]|2012l ; iBorriello et alj|2012j ). 

The most powerful testbed for these scenarios is the 
high-redshift Universe. Structure formation in WDM mod- 
els (or in any cosmological model with an equivalent power- 
spectrum c ut-off) is exponential l y sup pressed on small- 
scales (e.g., ISchneider etafl 12011 |201 j >. Since structures 
form hierarchically, these small halos are expected to host 
the first galaxies. If indeed dark matter were sufficiently 
"warm" , the high-redshift Universe would be empty. There- 
fore, the mere presence of a galaxy at high-redshift can set 
strong lower limits on the WDM particle mass. 

Due to their high luminosity, gamma-ray bursts (GRBs) 
constitute a remarkable tool to probe the high-z Universe 
and small-scale structures. Th ey provide a glimpse of th e 
first generations of stars (e.g.. Ide Souza et al.ll201ll . I2012T ) . 
as well as provide constraints on pri mordial non-Gau s sianit y 
|Maio et al.ll2012T ). As pointed out bv lMesinger et al.l (|2005h . 
the detection of a single GRB at z > 10 would provide very 
strong constraints on WDM models. 

Here we extend the work of iMesinger et all (|2005l ) by 
presenting robust lower limits on WDM particle masses, us- 
ing the latest Swift GRB data. The current data, includ- 
ing many redshift measurements, allows us to perform an 
improved statistical analysis by directly comparing the dis- 
tribution of bursts in various models as a function of red- 
shift. Furthermore, we make more conservative assump- 
tions throughout the analysis, such as normalizing the SFR- 
GRB ratio at high redshifts (thereby using a shorter, more 
accurate lever arm which minimizes modeling uncertainty), 
using an unbiased luminosity function and allowing the SFR- 
GRB ratio to evolve with redshift. Finally, we study the ef- 
fectiveness of future observations in improving the current 
constraints. 

Current limits on dark matter masses, m x , are mo- 
tivated by several observ ations. The Lym an-q forest im- 
plies m x > 1 keV (e.g. . IViel et alj 120081) and m vs > 8 



keV for sterile neutrinos (jSeliak et alj 2006; Bov arskv et al.l 
2009). Likewise, WDM models with a too warm candidate 
(m x < 0.75 keV) cannot simultaneously reprod uce the stel- 
lar m ass function and the Tully-Fisher relation (|Kang et al.l 
120131 ). Also, the fact t hat reionization occu rred at z > 6 im- 
plies m x > 0.75 keV (|Barkana et all 1200 ll ). However, all of 
these limits are strongly affected by a degeneracy between 
astrophysical (i.e. baryonic) processes and the dark matter 
mass. Our approach in this work is more robust, driven only 
by the shape of the redshift evolution of the z > 4 SFR. 
Furthermore, it is important to note that the SFR is ex- 
ponentially attenuated at high-redshifts in WDM models. 
Since astrophysical uncertainties are unable to mimic such 
a rapid suppression, probes at high-redshifts (such as GRBs 
and reionization) are powerful in constraining WDM cos- 
mologies. 

The outline of this paper is as follows. In §2 we discuss 
how we derive the dark matter halo mass function and SFR 
in WDM and CDM models. In §3 we derive the correspond- 



ing GRB redshift distribution. In §4 we discuss the adopted 
observed GRB sample. In §5 we present our analysis and 
main results. In §6, we discuss possible future constraints 
using a theoretical mock sample. Finally, in §7, we present 
our conclusions^ 



2 STRUCTURE FORMATION IN A WDM 
DOMINATED UNIVERSE 

Massive neutrinos from the standard model (SM) of par- 
ticle physics were one of the first dark matter candidates. 
However, structures formed in this paradigm are incompati- 
ble with observations. Other alternative dark matter candi- 
dates usually imply an extension of the SM. The DM par- 
ticle candidates span several order of magnitude in mass 
jBovarskv et al.ll2009h : axions with a mass of ~ 10~ 6 eV, 
first introduced to solve the problem of CP violation in 
particle physics, supersymmetric (SUSY) particles (grav- 
itinos, neutralinos, axinos) with mass in the range ~ eV- 
GeV, superheavy dark matter, also called Wimpzillas, (also 
considered as a possible solution to the problem of cosmic 
rays observed above the GZK cutoff), Q-balls, and sterile 
neutrinos with mass ~ keV range, just to c ite a few. For 
a revi ew about dark matter candidates see iBertone et al.l 
(2005). Two promising candi dates for warm dark mat- 
ter are the sterile neutrino f Dodelson fc Widrowl [1994: 



Shaposhnikov & TkachevU 2006) and gravitino is et al.l 
1984l;lMoroi et all 19931 ; iKawasaki et alll997i ; rPrimackll2003l ; 



Gorbunov et al.Tl2008l 



In WDM models the growth of density perturbations is 
suppressed on scales smaller than the free streaming length. 
The lighter the WDM particle, the larger the scale below 
which the power spectrum is suppressed. In addition to this 
power-spectrum cutoff, one must als o consider the r esidua l 
particle velocities. As described in iBarkana et al.l l|200ll ). 
these act as an effective pressure, slowing the early growth 
of perturbations. Bellow we describe how we include these 
two effects in our analysis (for more infor mation, please see 
IBarkana et"al]|200ll ; IMesinger et al.ll2005l ). 



2.1 Power spectrum cutoff 

The free-streaming scale, below whic h the linear pertur- 
bation amplitude is suppressed (e.g., IColombi et al.l [l996l : 
iBode et al.1l200ll ; IViel et al.ll2005l ). is given by the comoving 
scale 



Rr, 



0.11 



0.15 



1/3 



/ m x \ - 
VkeV/ 



-4/3. 



Mpc, 



(1) 



where il x is the total energy density contained in WDM par- 
ticles relative to the critical density, h is the Hubble constant 
in units of 100 kms _1 Mpc _1 , and m x is the WDM particle 
mass. 

The resulting modification of the matter power spec- 
trum can be computed by multiplying the CDM power 



1 Throughout the text, we use "conservative" to imply biasing 
the GRB distribution towards lower redshifts. This is conservative 
since it mimics the effects of WDM, thereby resulting in weaker 
constraints on the particle mass. 



2 Throughout t he paper, we adopt the cosmological parameters 
from WMAP9 llffinshaw et al.ll2012h . tt m = 0.264, C A = 0.736, 
n s = 0.97, erg = 0.8 and H = 71 km s" 1 Mpc _1 . 
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spectrum Pcdm ( fc) by an additional transfer function 
jBode et al.ll200ll l: 



PwDM(fc) = PcDM(fc) [1 + (ekf] 5M . 

where /i = 1.12 and 



(2) 



2.2 Effective pressure 

Structure formation in WDM models will be further sup- 
pressed by the residual velocity dispersion of the WDM par- 
ti cles, which delays the gro wth of perturba t ions. As shown 
in lBarkana et all |200ll ) and lMesinger et alj (120051 ). this "ef- 
fective pressure" is of comparable importance to the power 
spectrum cutoff in determining the WDM mass functions. 
The pressure can be incorporated in the halo mass function 
by raising the critical linear extrapolated overdensity thresh- 
old at collapse, 8 C (M, z). Using spherically symmetric hydro- 
dynamics simulations, and exploiting the a nalogy between 
the W DM effective pressure and gas pressure. iBarkana et al.l 
|200ll ) computed the modified 5 C (M, z). They showed that 
one can define an effective WDM Jeans mass, i.e., the mass 
scale below which collapse is significantly delayed by the 
pressure. We will denote this scale as: 



Mwdm * 1.8 x 10 1 



0.15 



1/2 



VlkeV/ 



(4) 



where we make the standard assumption of a fermionic spin 
1/2 particle. Note that equation (|4| diff ers from the original 
one proposed by IBarkana et al.l^200l] ) by a factor of 60. 
With this adjustment, we find that the collapsed fractions, 
F ca \\{z), simulated assuming a sharp mass cutoff at A/wdm 
are in good agreement with the full random- walk procedure 
(assum ing a more gradual rise in 5 C (M)) of IBarkana et all 
( 200ll ) (see Fig. [2J). The former approach has the advantage 
of being considerably faster and simpler. 



2.3 Collapsed fraction of haloes and cosmic star 
formation 

We use the Sheth -Tormen mass function, /st, 
jSheth fc Tormenl Il999h to estimate the number den- 
sity of dark matter haloes, ust(M, z), with mass greater 
than M: 



/st 



1 



ai(5f 



• exp 



(l\C 



2<T 2 



(5) 



where A = 0.3222, ai = 0.707, p = 0.3 and S c = 1.686. The 
mass function /st can be related to ust(M, z) by 

M dn ST (M,z) 



/st = 



p m diner -1 



(6) 



where p m is the total mass density of the background Uni- 
verse. The variance of the linearly extrapolated density field 
(j(M, z) is given by 



a(M,z) 



2tt 2 



r 

Jo 



k 2 P(k)W 2 (k,M)dk, 



(7) 



where b(z) is the growth factor of linear perturbations nor- 
malized to b = 1 at the present epoch, and W(k, M) is a real 
space top-hat filter. In order to calculate the CDM power 
spectrum Pcr>M(k), we use the CAMB cod<B For WDM 
models, we compute the power spectrum, Pwdm(&), with 
equation J2J. 

The fraction of mass inside collapsed halos, F co n(> 
-Mmm, 2), is then given by: 



F co n(z) = — 
p m 



dMMnsT(M,z) 



and the minimum mass is estimated as 

M min = Max[M gal (z),MwDM(m x )], 



(8) 



(9) 



where M ga i corresponds to the minimum halo mass capable 
of hosting star f orming galaxies. We u s e M ga i = M m i n from 
equation (11) in lSobacchi fc Mesingerl |2013af ). who present 
a physically-motivated expression for the evolution of M ga i 
which takes into account a gas cooling criterion as well as 
radiative feedback from an ionizingUV background (UVB) 
during inhomogeneous reionizatiorQ In other words, M ga i 
is set by astrophysics, whereas Mwdm(wx) is set by the 
particle properties of dark matter. In Fig. [T] we show these 
two limits as a function of redshift. As we shall see, our 
conclusions are derived from the regime where Mwdm > 
M ga i, and hence they are not sensitive to the exact value 
of M ga i. Another point worth highlighting is that models 
with m x > 3.5 keV have Mwdm(?ti x ) ^ M ga i up to z ~ 11. 
Thus, any constraint beyond ~ 3.5keV will no longer be 
robust, since galaxy formation can be suppressed even in 
CDM, bellow halo masses comparable to Mwdm- Therefore, 
observations at much higher redshifts are required to obtain 
tighter constraints. 

In figure we plot the fraction of the total mass col- 
lapsed into haloes of mass > M m i n , F co ii(> M m i n , z). The 
shaded region shows the collapsed fraction in CDM, with a 
range of low-mass cutoffs corresponding to virial tempera- 
tures 300 K < T v ir < 10 4 K. The other curves correspond 
to WDM particle masses of m x = 3.0, 2.5, 2.0, 1.5, and 
1. keV (top to b o ttom) . This figure is analogous to Fig. 2 
in lMesinger et all ((2005;) , serving to motivate equation Q. 
The fractions F co u(> M m i n , z), com puted according to th e 
full random walk procedure used in lMesinger et alT |2005). 
are shown as solid black curves, while the approximation of a 
sharp cutoff at Mwdm (eqf3| corresponds to the red dashed 
curves. 

The comoving star formation rate as a function of red- 
shift is assumed to be proportional to dF co \\/dt: 

dF co ii 



SFR(z) oc 



dt 



(10) 



The proportionality constant is irrelevant for our analysis, 
as it is subsumed in our normalization procedure below. 



3 http://camb.info/ 

4 Using M gal (2) from ISobacchi fc MesingeJ 1 1201331 1 for WDM 
models is not entirely self-consistent, since UVB feedback ef- 
fects should be smaller in WDM models. However, this effect is 
smaller than other astrophysical uncertainties. Most importantly, 
our conclusions are driven by models which at high redshift have 
Mwdm > Afgal , and are therefore insensitive to the actual choice 
of M gal (see Fig. 2). 
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Figure 1. The minimum halo masses capable of hosting star- 
forming galaxies. The so lid black line corresponds to the astro- 
physical limit, M ga i, from lSobacchi fc Mesingerl ll2013al) . The hor- 
izontal lines correspond to the cosmological cutoffs in WDM mod- 
els, Mwdm ( m x = 0.5,1,1.5,2,2.5,3,3.5 keV from top to bottom, 
respectively). 




Figure 2. Fraction of the total mass collapsed into haloes of 
mass > M m ; n as a function of redshift, F co n(> M m ; n ,z). The 
shaded region shows the collapsed fraction in CDM, with a range 
of low-mass cutoffs corresponding to virial temperatures 300 K 
< T v ; r < 10 4 K. The other curves correspond to WDM particle 
masses of m x = 3.0, 2.5, 2.0, 1.5, and 1.0 ke V (top to bottom ) . The 

~al] (2 



figure is an adapted version of Fig. 2 from lMesinger et al ] l |2005h 
(computed using their cosmology), serving to motivate equation 
Q. Values of F GO u(> M m j n ,z) computed acco rding to the full 
random walk procedure used in lMesinger et ahl l|2005h are shown 
as solid black curves, while the approximation of a sharp cutoff at 
Mwdm used in this work corresponds to the red dashed curves. 



3 THEORETICAL REDSHIFT DISTRIBUTION 
OF GRBS 

Under the hypothesis that the formation rate of long GRBs 
(LGRBs; dur a tion longer than 2 se c ) follows the SFR (e.g. , 
Totanil 1 19971; ICampisi et all |2010| ; iBromm fc Loebl |2006| ; 
de Souza et all l201ll ). the comoving rate of GRBs, ^grb, 



can be expressed as 

*grb(z) 



tn{l + z) a SFR(z) 



(11) 



Here £0 is a constant that includes the absolute conversion 
from the SFR to the GRB rate. The evolutionary trend 
described by a m ay arise from several mechanisms (e.g., 
iKistler et alll2009l ). with a possible explanation provided by 
the GRB preference for low-metallicity environment^. A 
metallicity threshold seems to provide a natural explanation 
for the observed value of a ~ 0.5 — 1 at low redshifts (z ^ 4; 
e.g. iRobertson fc Ellisl 120121 ). Such a metallicity threshold 
increases with redshift, whereas the characteristic halo mass 
decreases with redshift. In such a scenario, a value of a — 
would be appropriate for our high-redshift (z > 4) analysis. 
Nevertheless, we conservatively keep a as a free parameter. 

The observed number of GRBs in the range zi ^ z ^ zi , 
N{z\,z<z), can be expressed by 



N( Zl ,z 2 ) = K 



F 



2 rfiVGRB 

dz> 



dz , 



with 



dz 



*grb(z) 



At dV 

1 + z dz 



T{z) 



(12) 



(13) 



where the parameter K accounts for the efficiency of finding 
GRBs and measuring their redshift (e.g., area coverage, the 
survey flux limit, beaming factor of GRBs, etc) 0, dV/dz is 
the comoving volume element per unit redshift, At is the 



5 Since host galaxies of long duration GRBs are often observed 
to be metal poor, several studies have tried to connect the ori- 
gin of long GRBs with the metallicity of their progenitors (e.g., 
Meszaros!l2006l; IWoosley fc Bloomll2006l; [Salvaterra fc Chincarinil 
20071 : ISalvaterra et al.ll2009l ; ICampisi et al.ll201ll ). Such a con- 



ncction is physically-motivated since core-collapse models could 
not generate a long-G RB without the progenitor system hav- 
ing low metallicity (e.g., iHirschi et alj2005t lYoon fc Langerll2005l ; 
Woo slev fc Bloo m 2006). On the other hand, several authors re- 
port observations of GRBs in high metalli city environments (e.g., 
iLevesaue et al.ll20ld ; iKriihler et alj|2012h . suggesting that GRB 
hosts are not necessarily metal poor. Despite the apparent pref- 
erence of GRBs towards metal-poor hosts, there is no clear cutoff 
in metallicity, above which GRB formation should be suppressed. 
6 There are several selection effects known to mask the true GRB 
redshift distribution, e.g.: (i) the host galaxy dust extinction; (ii) 
the redshift desert (a redshift span, 1.4 < z < 2.5, in which 
it is difficult to measure absorption and emission spectra); (iii) 
Malmquist bias; and (iv) t he difference betw een redshift measure- 
ments techniques (e.g.. [Coward et al.ll2012h . We therefore expect 
K to be redshift dependent, and this evolution is subsumed in 
our parameter a above. Most of the above effects (e.g. obtaining 
GRB redshifts, dust extinction, Malmquist bias) result in biasing 
the observed sample towards low redshifts. Since we calibrate the 
proportionality between the GRB rates and SFRs at z ~ 3 — 4, 
we likely overestimate the efficiency in the redshift determination 
of z > 4 GRBs. Therefore, we expect that even a non-evolving K 
(i.e. our results for a = 0) would be a conservative assumption. 
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time interval in the observer rest fram^l, and X(z) is the 
integral over the GRB luminosity function (LF), p(L), 



l(z) = / p(L)dlogL. (14) 

ilogL Hm (z) 

To remove the dependence on proportionality constants, 
we construct cumulative distribution functions (CDFs) of 
GRBs over the redshift range zi < z < 2 max : 



N(< z|z max ) = 



N( Zi ,z) 



(15) 



N(zi, 2 max ) 

The expected number of observed GRBs in a given red- 
shift interval, for each combination 6 = {m x , a}, can be 
written as 



f Z2 At dV 

JV(*i,za;fl) = Co i » / (l + zrSFR(z;m x )-—-^l(z), 
J Z1 1 + z dz 



(16) 

where Co;fl is normalized so that each model recovers the ob- 
served GRB rate, N obB (z\, 22), at z ^ 4, 

iVob s (3,4) 



Co;9 = 74 



J^l + z^SFR(z;m^^l(z) 



(17) 



A' bs(3, 4) is equal to 24 (7) for SI (S2) samples respectively. 
Note that normalizing at a relatively high- z (z =3-4) is in- 
deed a conservative choice. If we had normalized at lower 
redshifts, say 2 ~ 1 — 2, the absolute number of GRBs pre- 
dicted by m x = 0.5 keV and CDM, between 2 ~ 3 — 4, 
would already diverge by a factor of ~ 2 for a values in 
the range — 2. In addition to being conservative, normal- 
izing at high-redshifts allows us to have a relatively short 
lever arm over which our simple scaling, SFR oc dF co n/dt, 
is presumed to be accurate. At lower redshifts, mergers and 
AGN feedback (missing from our model) are expected to be 
important in determining the SFR. Hence to normalize our 
CDFs, we chose the largest redshift at which our sample is 
reasonably large (see below). It is important to note however 
that our main results are based on the 2 > 4 CDFs, which 
unlike the predictions for the absolute numbers of bursts, do 
not depend on our choice of normalization. 



4 GRB SAMPLE 



Our LGRB data is taken from iRobertson fc Ellis! l|2012l) . 
cor responding t o a com p ilation from the sam p les presented 



Butler et all (120071); iPerlev et all (l2009h: [Butler et all 



ll2010h; ISakamoto et all (|201lD ; iGreiner et all (|201ll ), and 



iKriihler et all (|201ll ). It includes only GRBs before the end 
of the Second Swift BAT GRB Catalog, and is comprised of 
152 long GRBs with redshift measurements. It is important 
to consider the completeness of the sample. Several efforts 
have bee n made to construct a redshift-complete GR B sam- 
ple (e.g.. IGreiner et al|[20Tll ; ISalvaterra et alll2012l ). How- 
ever, to do so, many GRBs with measured redshifts are ex- 
cluded. Such requirements are even more severe for high-2 
bursts, which makes them of little use for our purposes. To 
explore the dependence of our results on a possible bias in 

7 Which in our case, corresponds to 5 yrs of observation time by 
Swift. Again, the value is not relevant for our purposes, since it 
is subsumed in the normalization. 
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Figure 3. Frequency (i.e. fraction in bin) of GRB luminosities 
for the z < 4 subsample used to construct the LF used for SI (see 
text for details). The red dashed line represents the best-fit LF. 



the GRB redshift distribution, we construct two samples: 
(SI) we use a luminosity function based on the 2 < 4 sub- 
sample (consisting of 136 GRBs); and (S2) we use a subsam- 
ple with isotropic-equivalent luminosities bright enough to 
be observable up to high redshifts (comprised of 38 bursts). 
The two samples are summarized in table [1] Since there is 
a degeneracy between a biased SFR-GRB relation and a 
redshift-dependent LF, we implicitly assume that any un- 
known bias will be subsumed in the value of the a parame- 
ter. 



4.1 Luminosity function sample (SI) 

The number of GRBs detectable by any given instrument 
depends on the specific flux sensitivity threshold and the 
intrinsic isotropic LF of the GRBs. In figure [3] we show the 
distribution of log-luminosities for 2 < 4 GRBs, which can 
be well described by a normal distribution 



P(L) 



p*e 



'(L-L'Y 



/2a 2 L . 



(18) 



logL iso /erg s 1 - L* 



log 10 5 



The values L 

1.06 and p* = 1.26 are estimated by maximum likelihood 
optimization. The luminosity threshold is given by 



Llin 



47TdiF Ull 



(19) 



where d_t is the lu minosi ty distance. Consistent with pre- 
vious works (e.g.[l]][2008), we set a bolometric energy flux 
limit Fnm = 1-15 x 10 _8 erg cm -2 s _1 for Swift by using the 
smallest luminosity of the sample. Due to Malmquist bias, 
our fitted LF is likely biased towards high luminosities. To 
minimize this problem and increase the sample complete- 
ness, we fit our LF using only 2^4 GRBs, which we show 
in Fig. |31 We reiterate that the Malmquist bias only serves 
to make our results even more conservative by predicting a 
flatter redshift distribution of GRBs. 
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Figure 4. Isotropic luminosities, Li so , of 152 Swift gamma- 
ray bursts as a function of z from the compilation of 
iRobertson fc EllidlipTj . The blue dot-dashed line approximates 
the effective Swift detection threshold (eg. 1191 . The black dashed 
horizontal line represents the luminosity limit of Li so > 1.34x 10 52 
ergs s — 1 , used to define our S2 subsample. 



4.2 Luminosity-limited sample (S2) 

Another approach, less dependent on the LF parametriza- 
tion and the Malmquist bias, is to construct a luminosity- 
limited subsample of the observed bursts bright enough to 
be seen at the highest redshift of interest. Assuming that 
the LF does not evolve with redshift, this subset would be 
proportional to the total number of bursts at any given red- 
shift. 

In figure [4] we show the redshifts and isotropic lumi- 
nosities of our entire sample. The dot-dashed blue line cor- 
responds to the effective Swift detection threshold. For our 
luminosity-limited sample, we only use GRBs with isotropic- 
equivalent luminosities Li so 1.34 x 10 52 ergs s _1 , which 
comprise all GRBs observable up to z ~ 9.4. Hereafter, all 
calculations will correspond to either the complete (LF de- 
rived) sample (SI), or the luminosity-limited sample (S2). 



5 OBSERVATIONAL CONSTRAINTS 

In this section, we test the WDM models by comparing the 
predicted absolute detection rates of bursts as well as the 
CDFs with the observed samples. We consider 3 different 
ranges of a: (i) a constant SFR-GRB relation, a = (case 
0, CO); (ii) —1 < a < 2 (case T CI); and (iii) a flat prior 
over — oo < a < oo (case 2, C2jf|. All cases are summarized 
in table [T] 



8 More precisely, C2 was run over the interval —3 < a < 12, 
which was more than sufficient to capture the likelihood decreas- 
ing to — > in the tails of the distribution (see Fig. [TJ) . 



Haiman, R. Perna, N. Yoshida 

Table 1. Set of cases considered in our analysis. 

Lli m (ergs s~ L ) a = — 1 < a < 2 — oo < a < oo 

L Hm 5 o sTco sTci sTc2 

L ]im > 1.34 x 10 52 S2C0 S2C1 S2C2 



5.1 Absolute detection rate of bursts 

In tables 1 2131 we present the absolute number of GRBs 
at high redshifts in CDM and WDM models with parti- 
cle masses of 0.5-3.5 keV, as well as the actual number in 
our sample observed with Swift. All models are normalized 
to yield the observed number of bursts at 3 < z < 4, as 
described in equation (|17p and the associated discussion. 

As expected, models with small WDM particle masses 
predict a rapidly decreasing GRB rate towards high red- 
shifts. This exponential suppression can in some cases be 
partially compensated by an increasing GRB-to-SFR rate 
(i.e. a > 0). For the S1C1 case, models with m x ^ 2.5 
keV show a good agreement with Swift observations for 
< a < 1, though values of a ~ 2 are a better fit to 
the observations at z > 8. For the case S2C1, all models 
with m x > 2.5 keV seem to be consistent with data for 
a ~ 1 — 2. In both cases, the two observed bursts in the in- 
terval 8 < z < 10, are already at odds with 1.5 keV < m x < 
2.5 keV models. Finally, we see that models with m x ^ 1 
keV predict a dearth of GRBs at z > 6, which is inconsistent 
with current observations. Extreme models with m x ~ 0.5 
keV already fail at intermediate redshifts (4 < z < 6), even 
for values of a as high as two. 



5.2 The redshift distribution of z > 4 bursts 

Although the absolute rate of bursts is the simplest predic- 
tion, it is dependent on the normalization factor between 
the SFR and GRB rate at 3 < z < 4. Hence, for the remain- 
der of the paper, we focus on comparing the theoretical and 
observed z > 4 CDFs. The CDFs are not dependent on nor- 
malization factors and are therefore more conservative and 
robust predictions. 

In figure[5]we plot the CDFs for CDM and WDM (under 
the assumption of a = 0), as well as the observed Swift 
distribution. The lighter the WDM particle, the sharper the 
CDF rise at low-z. There is a clear separation between CDM 
and WDM models with m x < 1.5 keV. Both the SI and S2 
samples (top and bottom panels respectively) show the same 
qualitative trends. 

As we saw above, the high- z suppression of structures in 
WDM models can be compensated for by allowing the GRB 
rate/SFR to increase towards higher redshifts. How degener- 
ate are these cosmological vs astrophysical effects? In figure 
[6] we show the CDF for m x = 0.5 keV for several values of 
a for S2 sample. The exponential suppression of DM halo 
abundances in this model is so strong, that an unrealistically 
high value of a ~ 15 is required to be roughly consistent with 
observations. Such a high value is ruled out by low-redshift 
observations, which imply a < 1 (e.g. | R.obertson fc Ellis! 
120121 ; iKistler et aT1l2009l ; iTrenti et alj|2012h . 
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Table 2. Absolute number of GRBs per rcdshift interval predicted by each model for S1C1 sample. 
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Table 3. Absolute number of GRBs per redshift bin predicted by each model for S2C1 sample. 



Model 




N(4,6) 






N(6,8) 






N(8,10) 






a = 


a = 1 


a = 2 


a = 


a = 1 


a = 2 


a = 


a= 1 


a = 2 


m x = 0.5 keV 


1.33 


1.65 


2.05 


0.01 


0.02 


0.03 


1.5 x 10~ 5 


3.1 x 10~ b 


6.6 x 10~ b 


m x = 1.0 keV 


4.10 


5.23 


6.70 


0.34 


0.59 


1.00 


0.01 


0.03 


0.06 


m x = 1.5 keV 


6.75 


8.78 


11.46 


1.27 


2.19 


3.76 


0.12 


0.26 


0.55 


m x = 2.0 keV 


6.47 


8.40 


10.94 


2.02 


3.52 


6.14 


0.31 


0.67 


1.44 


m x = 2.5 keV 


6.54 


8.49 


11.07 


1.86 


3.27 


5.77 


0.60 


1.30 


2.81 


m x = 3.0 keV 


6.57 


8.53 


11.12 


1.79 


3.13 


5.49 


0.74 


1.63 


3.58 


m x = 3.5 keV 


6.59 


8.56 


11.16 


1.80 


3.14 


5.51 


0.72 


1.58 


3.48 


CDM 


6.60 


8.56 


11.16 


1.80 


3.15 


5.51 


0.74 


1.63 


3.58 


Swift 


6 


6 


6 


3 


3 


3 


2 


2 


2 



5.3 Constraints from the redshift distribution of 
2 > 4 bursts 

To quantify how consistent are these CDFs with the ob- 
served distribution from Swift, we make use of two statistics: 
(i) the one-sample Kolmogorov-Smirnov (K-S) test; and (ii) 
a maximum likelihood estimation (MLE). Both tests are de- 
scribed in detail in appendix [Al 

The K-S test provides a simple estimate of the proba- 
bility the observed distribution was drawn from the under- 
lying theoretical one. We compute this probability, for fixed 
a first, for our models S1C0, S1C1, S2C0 and S2C1. Con- 
sistent with the more qualitative analysis from the previous 
section, models with m x < 1.0 keV are ruled out at 90% CL 
assuming -1 < a < 1. For a = (S1C0 and S2C0), the 
limits are even more restrictive and models with m x < 1.5 
keV are ruled out at 90% CL for both samples. 

So far, we have analyzed each model individually in or- 
der to quantify a lower limit on m x , given a single value 
of a. Using a \ 2 MLE (see appendix [X| allows us to com- 
pute posterior probabilities given conservative priors on a. 
Thus we are able to construct confidence limits in the two- 
dimensional, (m x , q) parameter space. The results for cases 
S1C2, S2C2 are shown in figure [7| at 68%, 95%, 99% CL. 
Both samples show the same qualitative trends, with the 
data preferring higher values of m x and CDM. Marginaliz- 
ing the likelihood over —3 ^ a ^ 12, with a flat prior, shows 
that models with m x ^ 1.6 — 1.8 keV are ruled out at 95% 
CL for S1C2 and S2C2 respectively. 



vations. We obtain constraints of m x > 1.6-1.8 keV. We 
now ask how much could these constraints could improve 
with a larger GRB sample, available from future missions? 
As a reference, we use the Sino-French space-based multi- 
band astronomical variable objects monitor (SVOMjf) mis- 
sion. The SVOM has been designed to optimize the synergy 
between space and ground instruments. It is forecast to ob- 
serve ~ 7 - 90 GRBs yr~* and ~ 2 - 6 GRB yr^ 1 at z ^ 6 
(see e.g., ISalvaterra et al.ll2008l '). 

We first construct a mock GRB dataset of 450 bursts 
with redshifts obtained by sampling the CDM, a — PDF 
given by equation (|13|) . This sample size represents an opti- 
mistic prediction for 5 yrs of SVOM observation^] (see e.g., 
ISalvaterra et al.l 120081 ). We then perform the MLE analysis 
detailed above on this mock dataset at z > 4. The resulting 
confidence limits are presented in Fig. [8] 

This figure shows that ~ 5 yrs of SVOM observations 
would be sufficient to rule out m x ^ 2.3 keV models (from 
our fiducial CDM, a = model) at 95% CL, when marginal- 
ized over a. This is a modest improvement over our current 
constraints using Swift observations. As already foreshad- 
owed by figures 1 and 2, as well as the associated discus- 
sion, it is increasingly difficult to push constraints beyond 
m x > 2 keV. On the other hand, the a constraint improves 
dramatically due to having enough high-z bursts to beat 
the Poisson errors. We caution that the relative narrowness 
around a — of the contours in Fig. [5] is also partially due 
to our choice of (CDM, a = 0) as the template for the mock 
observation. 



6 FUTURE CONSTRAINTS 

In the previous section, we have quantified the constraints 
on WDM particle masses using current Swift GRB obser- 



9 http://www.svom.fr 
The highest redshift in our mock sample is Zn 
the only event at z > 10. 



11.6, being 
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Figure 5. Cumulative number of GRBs for different values of m x 
compared with CDM predictions and Swift observations. The 
blue dotted line corresponds to m x = 0.5 keV, green dotted line 
to m x = 1.0 keV, red dotted line to m x = 1.5 keV, purple dotted 
line to m x = 2.3 keV, brown dotted line to m x = 2.5 keV, or- 
ange dotted line to m x = 3.0 keV, cyan dotted line to m x = 3.5 
keV, dark-green two-dashed line to CDM, black to the Swift ob- 
servations. Top Panel: sample S1C0; Bottom panel: sample 
S2C0. 



7 CONCLUSION 

Small-scale structures are strongly suppressed in WDM cos- 
mologies. WDM particle masses of m x ~ keV have been in- 
voked in order to interpret observations of local dwarf galax- 
ies and galactic cores. The high-redshift Universe is a pow- 
erful testbed for these cosmologies, since the mere presence 
of collapsed structures can set strong lower limits on m x . 
GRBs, being extremely bright and observable to well within 
the first billion years, are a promising tool for such studies. 

Here we model the collapsed fraction and cosmic SFR 
in CDM and WDM cosmologies, taking into account the 




4 5 6 7 8 9 



z 

Figure 6. Cumulative number of GRBs for m x = 0.5 keV as a 
function of the a parameter. Blue dotted line represents a = 0, 
green dotted line a = 3, red dotted line a = 6, purple dotted line 
a = 9, brown dotted line a = 12, orange dotted line a = 15, cyan 
dotted line a = 18, dark-green two-dashed line CDM and black 
line the Swift observations. 



effects of both free-streaming and effective pressure due to 
the residual velocity dispersion of WDM particles. Assuming 
that the GRB rate is proportional to the SFR, we interpret 
5 years of Swift observations in order to place constraints on 
m x . We conservatively account for astrophysical uncertainty 
by allowing the GRB rate/SFR to evolve with redshift as 
oc (1 + z) a . In order to fold completeness limits into our 
analysis, we used a low-z sample to estimate the intrinsic 
LF, or else restricted our analysis to a luminosity-limited 
subsample detectable at all redshifts. 

For each model (m x , a), we compute both the absolute 
detection rates and CDFs, at z > 4. A K-S test between the 
model and observed CDFs rules out m x < 1.5 (1.0) keV, 
assuming a = (< 2), at 90% CL. Using a maximum likeli- 
hood estimator, we are able to marginalize over a. Assuming 
a flat prior in a, we constrain m x > 1.6-1.8 keV at 95% CL. 
A future SVOM-like mission would tighten these constraints 
to m x > 2.3 keV. 

The strong and robust constraints we derive show that 
GRBs are a powerful probe of the early Universe. Their util- 
ity would be further enhanced with insights into their for- 
mation environments and their relation to the cosmic SFR. 
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(orange) and 99% (grey) probability. The asterisks correspond to 
the best-fit parameter combinations. Top panel: Sample S1C2; 
bottom panel: Sam ple S2C2 . Horizontal dotted lines represent 
the values a = lllshida et al.ll201ll; lElliott et alj|2012h . a = 0.5 
llRobertson fc Ellid2012ft . a = 1.2. jKistler et alj200Sh and a = 2 
for comparison. 
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APPENDIX A: STATISTICS 

Kolmogorov-Smirnov one-sample test 

A straightforward way to compare the data and WDM mod- 
els is to perform a one-sample Kolmogorov-Smirnov (K-S) 
test. The null hypothesis that the observed GRB redshifts 
are consistent with a model distribution can be evaluated by 
estimating a p-value, which corresponds to one minus the 
probability that the null hypothesis can be rejected. The 
K-S test consists in comparing the statistical parameter 



where F(z) and G(z) are the CDF for the theoretical and 
observed sample and sup is the the supremum of a totally or 
partially ordered set. We estimate the p-value via nonpara- 
metric bootstrap, which consists of running Monte Carlo 
realizations of the observed CDF using a random-selection- 
with-replacement procedure estimated from the data. This 
provides a histogram of the statistic D, from which a valid 
goodness-of-fit probability can be evaluated. The probabil- 
ity distribution function for each model is determined by 
equation (IB) . 

Maximum-likelihood method 

More formally, we can estimate the probability of parame- 
ters {a, m x } given the observed data using a Bayesian tech- 
nique. Assuming that our data is described by the prob- 
ability density function f(x;6), where a; is a variable and 
= {a, m x }. We want to estimate 0, assuming the data are 
independent. So the likelihood will be given by 



D = sup\F(z) - G{z)\, 



(Al) 



C(0\xi) oc Y[f(xi \B). 



(A2) 



Given the small number of observed bursts per red- 
shift biri J _A£_=_lJ) 2 we_use a Poisson error statistic^] (see 
e.g., ICampanelli et aJj|2012j for a similar procedure applied 
to galaxy cluster number count). Therefore, the likelihood 
function can be computed as 




where T, = N(zi, z i+1 ;0) and m = N obs (zi, z i+1 ). Thus, the 
\ 2 statistics can be written as 

X 2 (6) = -21n£, 



5 




(A4) 
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The value Az = 



1.5 is chosen to ensure at least 1 burst per 



bin. 
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